On Nonlinear Diffusion with Multiplicative Noise 
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Nonlinear diffusion is studied in the presence of multiplicative noise. The nonlinearity can be 
viewed as a "wall" limiting the motion of the diffusing field. A dynamic phase transition occurs 
when the system "unbinds" from the wall. Two different universality classes, corresponding to the 
cases of an "upper" and a "lower" wall, are identified and their critical properties are characterized. 
While the lower wall problem can be understood by applying the knowledge of linear diffusion with 
multiplicative noise, the upper wall problem exhibits an anomaly due to nontrivial dynamics in the 
vicinity of the wall. Broad power-law distribution is obtained throughout the bound phase. 
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PACS: 64.60.Ht, 02.50.-r, 47.20.Ky 

The diffusion equation with multiplicative noise has 
been used as a paradigm to describe a large class of 
nonequilibrium dynamical processes [y ranging from 
light propagation in random media to wealth fluctua- 
tions in economic systems. In this paper, we study the 
effect of nonlinearities, which inevitably occur in such 
processes due to the existence of characteristic scales of 
the diffusing fields. For instance, the growth of a popu- 
lation is bounded from above by a characteristic density 
determined by the competition of resources pi, and an 
investor's wealth can be taken as bounded from below by 
a fixed income source |^. At sufficiently coarse grained 
scales, these effects may be captured by the following 
nonlinear diffusion equation with multiplicative noise. 



dtn^W'^n-\-an-b n^+P + 7i(x, t) ri{x, t) 



(1) 



where a and b are constants, p specifies the degree of non- 
linearity, and ?7(x, t) is a Gaussian noise with zero mean 
and the variance (77(x,i)ry(x', i')> = "2.0 5{yi-yi')5{t-t'). 
Note that Eq. p) should be interpreted either in the Ito 
or Stratanovich sense H , depending on the details of the 
specific process being considered. We shall fix the non- 
linearity coefficient to be 6 = p from here on, so that for 
p > 0, the nonlinear term acts like a soft upper wall pre- 
venting n(x, t) from reaching large values even if a is large 
and positive, while for p < 0, the nonlinear term acts like 
a soft lower wall repelling n(x, t) from approaching zero 
even for large and negative a's. The specific choice of p 
depends on the symmetries and constraints of the system. 
The limits p — > ±oo mimic the effect of hard upper/lower 
walls and are also of interest. 

It will be convenient to cast the upper/lower wall prob- 
lems in a more "symmetric" form. By making use of the 
Cole-Hopf transform, h{x,t) = logn(x, i), we have 

dth = V^h + [Vhf + (a - ao) - pe^^ + 7/(x, t), (2) 

which is a Langevin equation similar to the KPZ equation 
describing the kinetic roughening of a growing interface 
h RH], with an additional drift term a — oq (oq = for 
Stratanovich dynamics and ao = D for Ito dynamics M), 



and an "exponential wall" at /i « 0. Thus a wall in 
the diffusion problem corresponds to a wall also in the 
interface problem. Eq. (||) shows that the sign of p de- 
termines the orientation of the wall while the magnitude 
of p describes the hardness of the wall. In the interface 
representation, it is clear that for large positive [nega- 
tive] a's, the system is pushed against the upper [lower] 
wall, while for large negative [positive] a's, the system is 
pushed away from the wall to h — — oo [h — +cxd], cor- 
responding to n = [n = oo]. A critical point separates 
the two regimes where the system "unbinds" from the 
wall. In this paper, we study the nature of this dynamic 
unbinding transition for walls characterized by different 
p's. We will show that while the soft and hard walls yield 
the same critical phenomena, the upper and lower walls 
yield two different universality classes. 

We start with a review of the exactly solvable "zero- 
dimensional" (single-site) version of Eqs. (13) and (0) [fj. 
Without spatial couplings, the stationary probability dis- 
tribution P{h) of Eq. (0) is simply 



P{h) oc exp 
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(3) 



where the second exponential merely enforces the con- 
straint of the wall, suppressing P for /i > [/i < 0] if 
p > [p < 0]. From (y), the stationary distribution of n 
is easily obtained using the relation P{n) — n^^ P(logn), 
yielding P(n) oc n^^"-"""'/^^-'^ e-x:p{-nP / D). Note that 
the distribution functions are normalizable as long as 
(a — ao)- > 0, and hence a power-law distribution in n 
is obtained for a large range of parameter values. [A 
similar result was obtained previously for the case of a 
hard lower wall [§U9|]; it was proposed as an explanation 
of the observed power-law distribution of wealth [Q .] 

As a -^ a J , the average h diverges according to Eq. (|^) 
as {h) ^ =p|a — ao|^^ for the upper/lower wall, indicat- 
ing the onset of an "unbinding transition" of h from the 
wall. This transition can be characterized quantitatively 
by singularities in moments of n. From P{n), we find for 
the upper wall problem (n™) ~ (a — ao)''" and (3m = 1 
for all integer m's. The distribution collapses towards 



P{n) = 5{n) as a ^ a^. For the lower wall, (n™) diverge 
as a — > a(^ since P(n) has no upper cutoff. It is conve- 
nient to characterize the phase transition in this case by 
monitoring n = 1/n. One finds (n™) ~ (oq — a)^", with 
Pm = 1 also. The symmetry between the upper and lower 
wall problems is evident from the h ^^ —h symmetry of 
Eq. (0) in the absence of spatial couplings. Note that the 
asymptotic critical properties given above are indepen- 
dent of the parameter p, indicating that in the vicinity of 
the unbinding transition, where (n) approaches zero or 
infinity (or as {h) -^ Too), the detailed form of the wall 
at finite n (or at /i w 0) is irrelevant. 

Time-dependent properties of the system can be ob- 
tained by solving the full Fokker-Planck equation [§,0. 
Qualitative features can be obtained alternatively by con- 
sidering the simpler Langevin equation for h{t). We illus- 
trate the solution by analyzing the problem with a lower 
wall {p < 0). For a < ao, h is confined to the range 
^ /i ^ 5h, where Sh ~ (oq — a)^^ is the scale of typical 
fluctuations in h. Right at the critical point a = oq, Sh 
diverges, i.e., 6h{t) ^ooast^oo. The form of Sh{t) is 
dictated by the equation of motion at the critical point, 
dth = -peP'' + ri{t). 5h{t) must be at leastoi 0{t^/'^) due 
to the random forcing ?/; the presence of the wall can only 
speed up the motion away from the wall. On the other 
hand, as h drifts far away from the wall, it should not be 
affected by the wall. Thus (5/i(i) ~ i^^'^ at the transition, 
with the distribution given by the scaling form 



P{h, t) = (Sh) ^ g [h/6h{t)] for ph < 0. 



(4) 



In (Q), the factor {Sh{t)) provides the proper normal- 
ization, and g{y) is a scaling function with the limiting 
behaviors g{y) « const for |y| <C 1 and g{y) — > for 
\y\ ^ 1. The distribution P{n,t) follows from (Q), yield- 
ing (ti™) '^ Sh ^ t^l"^ for all positive moments m. 

We now proceed to characterize the behavior of the 
system in arbitrary spatial dimension d. This is accom- 
plished by generalizing the above zero-dimensional pic- 
ture and utilizing the known properties of the KPZ equa- 
tion in finite dimensions [pi . The Cole-Hopf transforma- 
tion has already been exploited in Ref. [Q to derive some 
properties of the system with a soft upper wall (p > 0) . 
Some of the arguments there (but not all |12]) can be 
generalized also to the case of a lower wall. 

Close to the critical point, /i(x, i) is on average far from 
the wall, and the main source of fluctuation comes from 
the stochastic KPZ equation without walls (for all length 
scales below the correlation length ^). The scaling prop- 
erties of the KPZ equation can be described in terms of 
the dynamic exponent z alone g : the typical correlation 
time is T ~ ^^ and typical height fluctuation is Sh ^ ^~^ . 
These properties yield the scaling form for the additive 
renormalization {{'Vh)'^) of the bare drift (a — gq), with 



{{Vh)' 



0-1 — 02(^)1 where ai oc Z? is a constant and 



02(0 ~ (Sh/O'^ ~ ^-i/(22:-2)^ rpj^^g Q^g obtains the im- 
portant result ^ '^ |a — a^^'^ with Oc = ao — ai |13| ] 
and v — l/(2z — 2) for both the upper and lower wall. 



In particular, in one dimension where it is known ex- 
actly that z = 3/2, we have v = 1, while in d = 2, 
z ~ 1.60 yields i^ « 0.83. The exponents v and z can 
be used to relate the exponents /? [(3] and 0[9] which de- 
scribe the behavior of the "order parameter" n \n = n~^] 
close to the critical point for the problem with an upper 
[lower] wall. From the definitions (rt™) ~ [a — ad and 
(n™) ^ i"^™ at a = Oc, it follows that 6m — (3m/{i^z); 
simflarly, 6rn = P„Ji^^)- 

To find the exponents (3 and 9, we first perform a naive 
scaling analysis. Let us assume that the probability dis- 
tribution P{h, t) at the critical point is still given by 
the scaling form (Q) as in the zero-dimensional case, but 
with Sh{T) ^ t'^ where w = (2 - z)/z in d > 0. We 
then find (n™)uppoi. ^ ("■™}iowcr ~ Sh as before, yielding 
9„, = 9m =L0 = {2-z)/z and pm = ^,„ = (2-z)/(2z-2) 
for all to's. These exponents take on the value — 1/3, 
f3 = 1/2 in Id and 9 « 0.25 and f3 « 0.33 in 2d. 

In the above analysis, we have not distinguished be- 
tween systems with a lower or upper wall. On the other 
hand, it is evident that the equation of motion ^ for 
p > and p < are different in finite d due to the pres- 
ence of the nonlinear term (V/i)^. Thus the problems 
with upper and lower wall can no longer be mapped onto 
each other. In order to see whether the orientation of the 
wall is relevant to the critical behaviors, and whether the 
above scaling ansatz is correct, we perform a numerical 
simulation in one spatial dimension. We expose now the 
algorithm and the main results. 

To speed up the computations, we use a discrete model 
which, in the absence of walls, belongs to the KPZ univer- 
sality class . This model is similar to that of ballistic 
deposition: At each point a; of a one dimensional lattice of 
size L, we define a continuous height variable, ht{x). At 
every time step a new variable h'{x) = ht{x)-\-a + rjt{x) is 
obtained, with 7]t{x) being a random number uniformly 
distributed in the interval [0,1]. The value of h is up- 
dated simultaneously according to the rule ht+i[x) — 
max [ft,'(a; ± 1) -I- 7, ft.'(x, t)], where 7 = 0.1 is a fixed 
constant, and periodic boundary condition is applied. 
The parameter a controls the average drift of the sys- 
tem and is the only tunning parameter in this model. A 
hard wall a,t h — Q is introduced by including the ad- 
ditional rule ht{x) = min [/it(a;), 0] for an upper wall or 
ht{x) = max [/if (x), 0] for a lower wall. 

We describe first the results for an upper wall. The first 
step to studying critical properties is to locate the criti- 
cal point accurately. Starting with the initial condition, 
hi^{x) — for all x, we let the system evolve long enough 
so that a stationary state is reached {^ IQ^ time steps 
are typically required for system size oi L ^ l{fi). To 
compute any magnitude we average over up to 1000 inde- 
pendent runs. For a given size L, the critical point ac{L) 
is taken as the first value of a for which the steady state 
value (n) = (exp(/i)) becomes indistinguishable from zero 
as a is decreased. A plot of ac{L) vs 1/L is shown in 
Fig. 1(a). Linearity of the data suggests v ^ 1. Extrapo- 



lating to L — > oo, we get ac ~ —1.5750. The scaling of (n) 
and (n^) upon approaching the critical point ac{L) are 
shown in Fig. 1(b) for L = 400. We obtain the exponents 
/3i w /32 = 1.50 ± 0.15. In Fig. 1(b), we plot {n{t)) and 
(n2(i)) vs. i at a = a^ and find 6*1 « ^a = 1.10 ± 0.12. 
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FIG. 1. Results with an upper wall: (a) Size dependence 
of the critical point ac{L); (b) first two moments of the spa- 
tially-averaged n vs. a — ac in the steady state; (c) time 
dependence of (n) and (n^) at a = Oc; (d) {—h{t)) at a = Oc. 

The scaling of (h) at a = Qc is shown in Fig. 1(d). We 
find the exponent lo « 0.33. To compute the dynamic ex- 
ponent z, we perform a slightly different simulation p5[ : 
We take as initial condition a state in which h is a large 
negative constant at every lattice point except for the ori- 
gin where h = 0. The spreading of this "localized seed" 
is followed, and from the time evolution of the average 
size of the "infected region", we find z = 1.50 ± 0.05. All 
the measured exponents coincide, within numerical ac- 
curacy, with those computed previously 1 1 1 1 for Eq. (|l|) 
with p = 1,2, i.e., soft walls. This indicates that the 
introduction of an upper soft or a hard wall to a model 
belonging to the KPZ universality class yields the same 
critical behavior. As found in Ref. |11|, the exponents ly 
and z agree with the result of the scaling analysis. Also 
the scaling relation 9m — Pm/{vz) [ [ll|Jl5| is satisfied. 
However, the values of /3 and 9 are significantly different 
from those derived above using the naive scaling analysis. 

We next describe the results for a lower wall. Following 
the same analysis procedure for the "order parameter" 
(n) = (exp(— /i)), we obtain the critical point at ac ~ 
-1.5743, and 1/ « 1 (see Fig. 2(a)). From Fig. 2(b), we 
find /3i « /32 = 0.48 ± 0.05 @ and from Figs. 2(c) and 
2(d), we find Si « 6*2 = 0.32 ± 0.03, cj = 0.33 ± 0.04. The 
dynamical exponent z « 1.5 is again determined using 
the seed-spreading method (a suitable initial condition is 
h{x) = for all x except for the origin where h{Q) > 0). 
Note that in contrast to the case with an upper wall, 
all of the measured exponents here are consistent with 
the expected results based on the naive scaling analysis, 
including the values of (3 and 9. The clear difference 
between /?, 9 and /3, 9 indicate that problems with upper 
and lower walls belong to different universality classes. 
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FIG. 2. Results with a lower wall: (a) Size dependence 
of the critical point ac{L); (b) first two moments of the spa- 
tially-averaged n vs. flc — a in the steady state; (c) time 
dependence of (n) and (n^) at a = ac; (d) {h{t)) at a = ac- 

As the fluctuations in {h{t)) obey the same scaling law 
(uj « 0.33) for both the upper and lower wall problems 
(see Figs. 1(d) and 2(d)), differences in the scaling of 



{n"% 



and (n" 



/ lower 



must result from differences in 



the shape of the distribution P{h, t), or more specifically, 
in the form of the scaling functions glh/f^) defined in 
Eq. (H). Our numerical results suggest that the zero- 
dimensional behavior of g{y) « const for \y\ < 1 and 
g{y) — > for y ^ 1 is obtained in Id for the lower wall 
problem only. The numerically obtained forms of the 
scaling function g{h/t^^^) for the upper and lower wall 
are presented in Fig. 3(a) and Fig. 3(b) respectively. Note 
that for large \h\/t^^^, the distribution function drops 
off sharply (approximately exponentially) for both cases. 
Thus the assumed form for the scaling function g(jj) is 
satisfied for large y. However, the distributions are qual- 
itatively different for small values of \h\/t-^/^, reflecting 
qualitative differences in the interaction between h and 
the upper/lower wall at h = 0. 
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FIG. 3. Scaled probability distribution P{h, t) ■ f ^^^ for 10 
samples of size L = 3000, with (a) an upper wall [inset: log-log 
plot; straight line indicates {\h\/t^'^)'^\, and (b) a lower wall. 
The different curves are for t = 500, 2500, 5000, 10000. 

To appreciate the origin of this difference, let us con- 
sider the evolution of a piece of flat interface ft, — at 
the critical point a — Oc = oq ^ ai. The local drift rate 
(dth) « ((V/i)^) — oi is negative since h = and V/i = 
there. In the case of a lower wall, this negative drift 



has no effect since the interface cannot penetrate below 
the wall. A steady state results from a dynamic balance 
between the driving force {a — ac < 0) which "flattens" 
pieces of the interface against the wall, and the noise rj 
which roughens the interface. Thus, the situation is very 
similar to the zero-dimensional case, and we expect the 
exponent results obtained above from the naive scaling 
analysis to be valid for all d. 

The situation is very different for the case of an up- 
per wall, since the negative drift represents a significant 
force repelling the interface from the wall. The effective 
repulsion exerted by the upper wall is in fact long ranged: 
Suppose the interface is on average at a distance h from 
the wall. Then interfacial roughness is "interrupted" by 
the wall at a scale ^ ~ /i^ in Id, leading to an effective 
repulsive force {dth) ~ — l/f"^ ~ —l/Ji^. Thus, a steady 
state can only be achieved by keeping the interface away 
from the wall. This necessarily leads to the depletion 
of P{h,t) for small \h\. In analogy to the problem of a 
random walker in the vicinity of an absorbing wall, we 
expect P to have a power law form, P{h,t) ^ {\h\/t^'^)"' 
for \h\ <^ t^/'^. Our numerical results (Fig. 3(a) inset) are 
consistent with this power law form, with the exponent 
(T ~ 2. The new form of P for the upper wall allows 
us to obtain moments of n in terms of the exponent a. 
We find 6m = {1 + o) ■ lo w 1. Our numerical result 
9m ~ 1.1 ± 0.1 (Fig. 1(c)) is consistent with this expo- 
nent relation. The anomalous exponents can in principle 
be computed by analyzing the effect of nonlinearity along 
the line of Ref. |ll[ , but starting with the exact knowl- 
edge of the KPZ equation in d = 1. 

Away from the critical point, the form of P{h) is 
simply obtained by replacing 5h{t) ^ f^ with \a — 



-v{2-z) 



-1/2 



1. Our numerics 



for the lower wall problem (Fig. 3(b)) suggests that 
P{h) oc exp [— cft,(ac — a)^'^] where c is a nonuniver- 
sal constant. This exponential distribution of h (which 
we have also checked directly by numerics) leads to a 
power law distribution of n as in the zero-dimensional 
case. Again, we expect this result to be valid for all 
d. For the upper wall problem, the depletion^effect de- 
scribed by the extra factor [\h\{a — acY^"^]" in P{h) leads 
to logarithmic corrections to the power law form of P{n). 
In summary, nonlinear diffusion with multiplicative 
noise has been analyzed in terms of the dynamic unbind- 
ing of an interface from a wall. Two distinct universality 
classes are obtained for the upper and lower walls in fi- 
nite spatial dimensions. The upper wall problem, which 
describes various saturation effects, exhibits anomalous 
scaling due to nontrivial interaction with the wall. In 
contrast, the lower wall problem, which models the in- 
teraction of a growing interface with a substrate and the 
fluctuation of wealth in economic systems, can be under- 
stood by combining and applying the knowledge of the 0- 
d problem and the KPZ equation. It is interesting to note 
that the lower wall problem, in particular the numerical 
algorithm used here for the hard lower wall, resembles 



the algorithms used for the "local, gapped" alignment 
of DNA sequences |l^,|l^. (Similarly, the 0-d problem 
is analogous to gapless local alignment [1^].) Divergent 
critical fluctuations described in this work have a signifi- 
cant effect on the optimization of sequence alignment and 
will be discussed in detail elsewhere. 
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